Libraries


library(terra)
library(tidyterra)
library(ggplot2)
library(dplyr)
library(scales)

# Get the data
library(geodata)

Get the data


# I use to save this data on my local machine to speed up the code

mycachedir <- "~/R/mapslib/misc"
r <- elevation_30s("COL", path = mycachedir)

# For better handling we set here the names
names(r) <- "alt"

# Quick look
autoplot(r) +
  theme_minimal()

Hillshading

See also https://dominicroye.github.io/en/2022/hillshade-effects/


## Create hillshade effect

slope <- terrain(r, "slope", unit = "radians")
aspect <- terrain(r, "aspect", unit = "radians")
hill <- shade(slope, aspect, 45, 270)

# normalize names
names(hill) <- "shades"

# Hillshading, but we need a palette
pal_greys <- col_greys <- grey(1:100 / 100)

ggplot() +
  geom_spatraster(data = hill) +
  scale_fill_gradientn(colors = pal_greys, na.value = NA)

For blending we follow a different approach


# Use a vector of colors


index <- hill %>%
  mutate(index_col = rescale(shades, to = c(1, length(pal_greys)))) %>%
  mutate(index_col = round(index_col)) %>%
  pull(index_col)


# Get cols
vector_cols <- pal_greys[index]

# Need to avoid resampling
# and dont use aes

hill <- ggplot() +
  geom_spatraster(data = r, fill = vector_cols, maxcell = Inf)

hill

A note on colors


# Try some options, but we need to be aware of the values of our raster

r_limits <- minmax(r) %>%
  # Add some extra values to max
  as.vector() %>%
  # Minimum on 0
  pmax(0)

# Rounded to closest 500
r_limits <- ceiling(r_limits / 500) * 500

# Compare
minmax(r) %>% as.vector()
#> [1]  -80 5452
r_limits
#> [1]    0 5500


# Now lets have some fun with scales from tidyterra

elevt_test <- ggplot() +
  geom_spatraster(data = r)

# Create a helper function

plot_pal_test <- function(pal) {
  elevt_test +
    scale_fill_hypso_tint_c(
      limits = r_limits,
      palette = pal
    ) +
    ggtitle(pal) +
    theme_minimal()
}

plot_pal_test("etopo1_hypso")
plot_pal_test("dem_poster")
plot_pal_test("colombia_hypso")
plot_pal_test("pakistan")
plot_pal_test("utah_1")
plot_pal_test("wiki-2.0_hypso")